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ABSTRACT 

Multipath propagation in the interstellar medium distorts radio pulses, an effect predominant for 
distant pulsars observed at low frequencies. Typically, broadened pulses are analyzed to determine the 
amount of propagation-induced pulse broadening, but with little interest in determining the undistorted 
pulse shapes. In this paper we develop and apply a method that recovers both the intrinsic pulse shape 
and the pulse broadening function that describes the scattering of an impulse. The method resembles 
the CLEAN algorithm used in synthesis imaging applications, although we search for the best pulse 
broadening function, and perform a true deconvolution to recover intrinsic pulse structre. As figures of 
merit to optimize the deconvolution, we use the positivity and symmetry of the deconvolved result along 
with the mean square residual and the number of points below a given threshold. Our method makes 
no prior assumptions about the intrinsic pulse shape and can be used for a range of scattering functions 
for the interstellar medium. It can therefore be applied to a wider variety of measured pulse shapes and 
degrees of scattering than the previous approaches. We apply the technique to both simulated data and 
data from Arecibo observations. 

Subject headings: Pulsars, Interstellar Scattering, Data Analysis, CLEAN Technique 



1. INTRODUCTION 

Multipath propagation in the interstellar medium (ISM) 
at radio wavelengths is caused by diffraction from irregu- 
larities in electron density on scales smaller than the Fres- 
nel scale, VXD ~ 10 11 cm (where A is the wavelength of 
observation and D is the extent of the scattering medium) . 
One consequence is that pulses from radio pulsars are dis- 
torted by the differences in arrival time of different por- 
tions of the scattered wavefronts. The resultant pulse 
broadening (e.g. Williamson 1972) is usually quantified in 
terms of the pulse-broadening time, r<j, which is related 
to the scattering measure SM, where SM is the integrated 
value of C„, the coefficient of the wavenumber spectrum of 
electron density irregularities. Pulse broadening is impor- 
tant information for quantifying the Galactic structure of 
electron density fluctuations. But it is also a hindrance in 
studying the intrinsic properties of pulsars that are mani- 
fested in the undistorted pulse shape. 

For a Kolmogorov wavenumber spectrum with small in- 
ner scale (e.g. Rickett 1990; Cordes & Lazio 1991; Arm- 
strong, Rickett & Spangler 1995), the pulse broadening 
time scales as 

r d « 1.10 ms W T SM 6/5 iy- 22 / 5 D, (1) 

where v is the frequency of observation (GHz), D is the 
pulsar distance (kpc), SM has units kpc m -20 / 3 , and W T 
is a geometric factor that accounts for the distribution 
of C 2 along the line of sight (e.g. Cordes & Rickett 1998). 
Pulsars with larger dispersion measures (DM) tend to have 
larger scattering measures. The effect is thus predominant 
in observations of distant pulsars and at lower frequencies, 
and it limits the utility of such scattered pulse profiles for 
the study of pulsar magnetospheres. 

The observed pulse y(t) can be modeled as a convolu- 
tion of the intrinsic pulse x(t) with the pulse-broadening 
function (PBF) g(t) and the net instrumental resolution 



function, r(t): 

y(t) = x(t)®g(t)®r(t); (2) 

the PBF g(t) is the response of the interstellar medium 
(ISM) to a delta function. The exact form of the PBF 
depends on the spatial distribution of scattering material 
along the line of sight and on its wavenumber spectrum 
(e.g. Williamson 1972, 1973; Cordes & Rickett 1998; Lam- 
bert & Rickett 1999; Cordes & Lazio 2001; Boldyrev & 
Gwinn 2002). 

In general, reconstruction of the original pulse shape 
from the scattered data is complicated by the fact that 
the PBF g(t) is unknown. However, classes of possible 
PBFs corresponding to different distributions of scattering 
media are known. 

Independent determination of the PBF requires knowl- 
edge of the wavenumber spectrum and distribution of scat- 
tering material along the line of sight. Alternatively, a 
priori knowledge of the intrinsic pulse shape could allow 
determination of the PBF. In general, however, both x(t) 
and g(t) are unknown; hence the utility of scattered pulse 
profiles is limited to measuring t c i through a template fit- 
ting method, which assumes a form for the intrinsic pulse 
shape and a model for the ISM, rather than determining 
these quantities from the data. The assumed shapes can 
be incorrect, leading to estimates for which are only 
approximate, and have poorly quantified uncertainties. 

In this paper, we describe a method for reconstruction 
of the intrinsic pulse shape as well as determination of 
the PBF. The basic idea is derived from the CLEAN algo- 
rithm extensively used in image reconstruction of synthesis 
imaging data (Hogbom 1974; Schwarz 1978). Our one- 
dimensional application is similar to the usage of CLEAN 
in spectral analysis (Roberts, Lehar & Dreher 1987). In 
§2, we outline the method and its practical realization, 
and compare it to other methods that have been used for 
this application. This CLEAN-based method is demon- 
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strated using simulated data as well as some recent data 
from Arecibo observations, and the results are discussed 
in §3. Finally, we present our conclusions in §4. 

2. METHOD AND PRACTICAL REALIZATION 

Before outlining the CLEAN-based approach, we briefly 
describe other methods that have been used to estimate 
the pulse-broadening time from scattered pulse profiles. 
In general, a model is fitted to the observed pulse profile, 
y (t) , where the model, y m (t) , is the convolution of a 
template profile, representing an estimate for the intrinsic 
pulse shape, with an impulse response, g(t), that describes 
pulse broadening from the ISM for the particular line of 
sight. 

2.1. Conventional Frequency- Extrapolation Approach 

Typically, the intrinsic pulse shape at the frequency of 
interest is extrapolated from high-frequency observations 
where the pulse broadening is assumed to be negligible. In 
addition, the PBF shape is also assumed known, often as 
a one-sided exponential function, 

g(t) =T d - 1 cxp(-t/T d )U(t), (3) 

where U (t) is the unit step function. As discussed below, 
this form applies to a thin scattering screen for which the 
density irregularities follow a square-law structure func- 
tion. The pulse broadening time T d is estimated by mini- 
mizing x 2 ( as standardly defined) as a function of r<j, the 
sole parameter in this problem. Results from the method 
thus rely on use of an empirical pulse shape derived from 
radio frequencies quite different from that at which the 
scattering is prominent. However, it is well known that 
profile evolution with frequency can be quite complex for 
some pulsars (e.g. Hankins & Rickett 1986), and extrapo- 
lation to the frequency of interest will be in error. It is also 
well known that the PBF is a one-sided exponential only 
for a medium described as a thin screen with a square-law 
structure function (e.g. Cordes & Rickett 1998). While 
thin screens may exist along some lines of sight, a more 
reasonable default model is one where the medium is Kol- 
mogorov (leading to a non-square law structure function) 
and fills a significant fraction of the line of sight. Despite 
such shortcomings, the method has been used to obtain 
estimates of T d for a large number of pulsars (e.g. Lohmer 
et al. 2001). 

Variations on this approach include use of non-exponential 
pulse broadening functions and detailed estimation of the 
frequency dependence of the intrinsic pulse shape, for ex- 
ample by modeling individual components of the pulse 
shape and scaling each with frequency. The essence of the 
method is the same, however, in that measurements at 
other frequencies are used to estimate the intrinsic pulse 
shape at the frequency of interest. 

2.2. Fourier Inversion 

From Eq. 2, the Fourier transform of y(t) is related to 
that of x(t) as Y(f) = X(f)G(f)R(f), where X, G and R 
are the respective Fourier transforms of the factors on the 
right-hand side of Eq. 2. Adopting a particular form for 
g(t) and knowing r(t), one can, in principle, perform the 
deconvolution by calculating X(f) = Y(f)/G(f)R(f) and 
inverse-transforming to obtain the intrinsic pulse shape. 



This method was used by Weisberg et al. (1990) and 
Kuz'min & Izvekova (1993) with a one-sided exponential 
PBF. Kuz'min & Izvekova pointed out that the results de- 
pend on the assumed form of the PBF but did not present 
any results using alternative forms. 

2.3. A CLEAN-based approach 

In this paper, we propose an alternative method to real- 
ize both objectives of recovering the intrinsic pulse shape 
and determining the shape and characteristic time scale 
of the broadening function. The method is based on the 
CLEAN algorithm, an iterative method for deconvolving 
the instrumental point source response function from an 
intcrferometric image (Hogbom 1974; Schwarz 1978). In 
the present case, the problem is much simplified as the 
data are one-dimensional. The method involves inverting 
the convolution in Eq. 2. By analogy with imaging appli- 
cations, the measured, scattered pulse is equivalent to the 
"dirty map" ; and the combination g(t)®r(t) (the PBF con- 
volved with the instrumental response) is equivalent to the 
"dirty beam" . In synthesis imaging applications, CLEAN- 
ing involves the iterative subtraction of scaled copies of the 
dirty beam from the dirty map at locations corresponding 
to real features until the residual map is indistinguish- 
able from noise, followed by restoration of the subtracted 
CLEAN components. Here, we adapt the method to (one- 
dimensional) scattered pulse profiles. 

While this method is inspired by CLEAN as used in syn- 
thesis imaging, we emphasize that there are some impor- 
tant differences between the two algorithms. Specifically, 
in imaging applications, the dirty beam is assumed to be 
perfectly known, and sets an upper limit on the achiev- 
able resolution. In the algorithm outlined here, we are 
performing a true deconvolution, and the resolution of the 
deconvolved pulse is determined by the response function 
of the instrument, rather than the width of the PBF. Ad- 
ditionally, the PBF is not known a priori, and we perform 
a search in both parameter space and function space for 
the best pulse-broadening function, as outlined below. 

We consider the measured pulse to consist of delta func- 
tions, each convolved with the PBF, g(t), and the resolu- 
tion function, r(t). In reality, radiometer noise adds to 
the measured quantity and the intrinsic pulse shape is it- 
self perturbed by pulse-to-pulse fluctuations intrinsic to 
the pulsar. The CLEAN process identifies the amplitude 
and location of a CLEAN component (CC) by finding a 
suitable maximum in the residual profile in each iteration. 
Initially, the residual profile is the measured profile, y(t). 
The amplitude of the CC is the identified maximum mul- 
tiplied by a loop gain, 7 . The CC y c (t) is convolved with 
g(t) and r(t) and subtracted from the measured pulse pro- 
file to yield the residual pulse profile, Ay(t): 

Ay(U ) = y(U ) - [y c (t) ® (g(t) ® r(*))] t=tj , i = l, N, 

(4) 

where y c (U ) = 7 {max[y(t)]}<5(t — t ), and N is the num- 
ber of pulse phase bins. Ay(t) for a given iteration be- 
comes y(t) for the next iteration, and the process is re- 
peated until the residual subtracted pulse, Ay(t), is com- 
parable to the off-pulse rms noise level 1 . Upon termina- 
1 

The termination criterion is not well-defined for the CLEAN tech- 
nique. 
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tion of the iteration, n c CLEAN components are identi- 
fied: Cj,tj, j = 1, . . . , n c denote the amplitudes and times 
of these CLEAN components (i.e., the collection of delta 
functions at the end of the CLEAN process). 

The restored pulse shape is built from the ensemble of 
CLEAN components, Cj, convolved with a 'restoring func- 
tion' p{t) , details of which we will discuss later in this sec- 
tion: 



g(t) = (7r 5 r//8t 5 ) 1/2 exp(-7r 2 T d /4i) 



(7) 



y r (t) = J^C j 5(t-t j )®p{t). 



(5) 



At the end of the CLEAN process the residual noise, 
Ay(t), is added to the convolved CCs to obtain the final 
restored pulse profile. 

Examples of the CLEAN process are shown in Figure 1 
for simulated data (left-hand panels) and real data (right- 
hand panels). The figure shows the measured pulse shapes, 
the chosen PBF (a one-sided exponential as in Eq. 3), 
the CCs, and the deconvolved profile. Compared to the 
traditional least-squares or frequency-extrapolation meth- 
ods (§2.1), the approach to the pulse-broadening problem 
outlined here is much easier to apply because it is non- 
parameteric. Further, it does not make any assumption 
about the intrinsic pulse shape or its evolution with fre- 
quency, and thus can potentially yield useful information 
about the intrinsic pulse shape. 

Below we elaborate on some aspects concerning the 
practical realization of the method. 

2.3.1. Loop Gain 

Small values of the loop gain (7) are employed in imag- 
ing applications. The performance of the method is ex- 
pected to change only marginally for sufficiently small 
values (say, 7 <C 0.1). Larger values tend to cause 
over-subtraction and consequent lack of convergence. We 
adopted 7=0.05 as a reasonable choice for our applica- 
tion. Thus, in each iteration, a delta function with an 
amplitude of 5% of the profile peak is convolved with the 
chosen PBF and the instrumental response function r(t), 
and then subtracted from the peak of the residual pulse 
profile. 

2.3.2. Pulse Broadening Functions 

Unlike radio synthesis imaging applications, where the 
dirty beam is often known exactly as the response of the 
instrument to a point source, the exact functional form for 
the PBF of the ISM is not known. It depends on both the 
distribution and the wavenumber spectrum of scattering 
material along the line-of-sight to the pulsar, neither of 
which is known in detail. The PBF is particularly sensi- 
tive to the spatial distribution of C 2 . For instance, in the 
simplest case of a thin slab of scattering material located 
in between the pulsar and the Earth, the PBF is character- 
ized by a sharp rise and exponential fall-off (cq. 3), whereas 
the PBFs of more realistic models (such as a thick slab or a 
uniform medium) are characterized by rounded edges and 
slower fall-offs, with different rise and peak times. For me- 
dia with a square-law structure function (Cordes & Rickett 
1998; Lambert & Rickett 1999) the PBF for a thin screen 
is given by Eq. 3 and for thick and uniform media the 
PBFs are given by Williamson (1972, 1973): 



g(t) 



(7TT d /4t 3 ) 1/2 CX P (~Tr 2 T d /m) 



(6) 



where the time constant r d is the pulse-broadening time. 

For other media, such as those with a Kolmogorov 
wavenumber spectrum having different inner scales (Lam- 
bert & Rickett 1999) or media with non-Gaussian density 
fluctuations (Boldyrev & Gwinn 2002) (but still having 
uniform statistics transverse to the line of sight), the PBFs 
are qualitatively similar to those for media with square- 
law structure functions in that thin screens show sharp 
rise times while extended media have slower rise times. 

Nonuniformities of the scattering transverse to the line 
of sight, particularly any truncation of the scattering 
screen, produce PBFs truncated for times t > t max , where 
imax is the arrival time corresponding to the edge of the 
screen (Cordes & Lazio 2001). As an example, a circular 
disk centered on the line of sight to a source will produce 
a PBF of the form 

"" = r^(- ) P "---"- <8) 

where t max = Cu, C = Cx/ 2<T 8. where max and a g 
are the maximum and RMS scattering angles, respectively. 
Cordes & Lazio (2001) also discuss PBFs for several other 
configurations. For example, in the case of scattering from 
a long, narrow filament, # max will be much larger than ag , 
and consequently the scattered image becomes elongated, 
and the PBF tends toward the form cx t -1 / 2 exp(— t). For 
a large screen with elongated irregularities, the PBFs will 
be strongly frequency dependent, and will be characterized 
by multiple peaks that align at different frequencies. 

In applying our algorithm, one can perform a grid search 
over various forms of the PBF as well as the relevant pa- 
rameters and thus determine (a) the appropriate model 
of the PBF that best describe observations, and (b) the 
best fit Td that yields meaningful intrinsic pulse shapes. 
The figures of merit for determining the best fit case is 
discussed later in this section. 

2.3.3. Instrumental Response 

The response function due to instrument can be de- 
scribed as the convolution of the functions that charac- 
terize various instrumental effects relevant in determining 
the effective pulse phase resolution, including (i) disper- 
sion smearing (r dm ) , (h) profile binning (r p t ) , (iii) the 
backend time averaging (r av ), and (iv) additional post- 
detection time averaging (r pd ) . The response function of 
the instrument, or the effective resolution function of the 
data, is therefore given by 

r(t) = r dm (t) <g> r pb (t) <g> r av (t) <g> r pd (t). (9) 
The contributing factors in the above equation, in partic- 
ular those which describe effects (i) to (iii), can be treated 
as approximately rectangular. Hence the response func- 
tion will be approximately "trapezoidal" for cases where 
one or two factor(s) dominate the effective time resolu- 
tion of the data. It will be more "Gaussian-like" due to 
the convolution of similar-sized rectangle functions when 
there arc multiple contributions that determine the time 
resolution. We form the resolution function according to 
eq. 9, where the convolution is performed with a time res- 
olution that is much larger than the narrowest factor; the 
resultant function is then re-sampled to the resolution of 
the actual data. 
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Fig. 1. — Sample plots that illustrate the CLEAN process for deconvolving the interstellar pulse broadening. Left-hand side: simulated 
data. Right-hand side: data for PSR J1852+0031 at 1475 MHz. (a) scattered pulse shapes (b) model PBFs (c) CLEAN components (d) final 
restored (CLEANcd) pulse shapes. 



2.3.4. Restoring function 

The CLEANed pulse shape is obtained by convolving 
CCs with a restoring function, pit), equivalent to the 
CLEAN beam in imaging applications. The choice of a 
restoring beam is driven by several ideas (see Cornwell, 
Braun & Briggs (1999) for a detailed discussion). Follow- 
ing the analogy of imaging applications, where the CLEAN 
beam is approximated as a Gaussian fitted to the central 
part of the dirty beam, we adopt an equivalent Gaussian, 
r g (t), whose width reflects the effective resolution of the 
actual data (given by the width of r(t) as discussed in 
the previous section). The choice of restoring function, 
p(t) = r g (t), de-emphasizes the higher frequency com- 
ponents that are sometimes spuriously generated by the 
CLEAN algorithm when directly using the response func- 
tion r(t) as the restoring function. The restoring function 
is normalized to unit area in order to ensure conservation 
of flux in the inversion process. 

2.3.5. Figures of Merit 

We now discuss the figures of merit for determining the 
best fit case, i.e., when the observed pulse profile is op- 
timally compensated for interstellar broadening. Assume 
that we know the shape of the PBF but not its character- 
istic scale, Td- If we choose Td too large, we will tend to 
oversubtract pulsed flux and thus introduce negative go- 
ing features in the deconvolved profile. Conversely, if too- 
small a value for Td is used, the deconvolved profile will 
display residual asymmetry due to un-deconvolved scat- 
tering. Therefore two figures of merit that are appropri- 



ate are: (1) positivity and (2) minimum asymmetry of 
the resultant pulse. In practice we numerically require 
the residual noise (i.e., measured— modeled pulse) to show 
zero mean and the same rms everywhere. We use all these 
criteria to choose the correct PBF and characteristic time 
scale, Td- 

We define a parameter (f r ) that measures positivity as 

N 

•^ = 7T^E[ A ^)] 2C/ Ay (10) 
JVCT °ff' i=1 

where 

U A y = U(-Ay(ti)-xa oS ), (11) 
Ay(ti ) is the residual noise at the end of the CLEAN pro- 
cess, and m denotes a 'weight' of order unity. The unit 
step function U turns on when the residual Ay is signifi- 
cantly below the off-pulse noise level (i.e. when oversub- 
raction is caused by using a PBF that is too broad). We 
adopted x — 3/2 for our analysis so that distortions of the 
pulse more negative than 3cr/2 yield a penalty manifested 
as a larger value of f r . To minimize profile asymmetry, we 
minimize the skewness parameter of the CLEANed pro- 
file, y r (t), by explicitly minimizing the skewness (r) of the 
CCs: 



where 

n c 

n = ^ > (13) 

i=i 
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t = ^h — . a*) 

i=i 

and Ci,ti, i = 1, . . . ,n c denote the amplitudes and times 
of the CLEAN components, as defined before. 

In addition to the symmetry and positivity constraints, 
we also maximize the number of residual points in the 
on-pulse window that were consistent with the noise in 
an off-pulse window. Specifically we count the number of 
points Nf that satisfy \yi — (y s) | < 3cr ff , where (y g) and 
(7 ff represent the mean and rms of the off-pulse region. 
Thus the figures of merit we employ include (a) maximum 
Nf, (b) minimum residual noise, (c) minimal skewness, 
and (d) minimum positivity parameter f r . 

Sample plots illustrating the figures of merit are shown 
in Figure 2 for the simulated and actual data of Figure 1. 
The best fit model will correspond to minimal values of 
both f r and T, or equivalently to a minimum of the com- 
bined parameter, f c = (T + f r )/2, as well as conforming to 
a maximum of Nf and a minimum of the residual noise rms 
(c ffc ) of the CLEANed pulse. Economy in the number 
of CCs required to fit the observed pulse profile may also 
provide a means to discriminate between different forms of 
the PBF (and hence, different distributions of scattering 
material). 

As the deconvolution process removes smearing due to 
interstellar scattering, the resultant CLEANed pulse is 
narrower and the peak signal to noise (S/N) is higher, 
the increase being determined largely by the ratio of Td to 
pulse width. This is evident from examples of simulated 
and real data shown in Figs. 1 and 3-6; examination of 
Figure 1 suggests that the sum of CCs matches closely the 
area of the scattered pulse. 

3. FURTHER EXAMPLES 

We illustrate the application of the above-described 
method using additional simulated and real data. In doing 
so, we show the sensitivity of the method for finding both 
the shape of the pulse broadening function and the appro- 
priate Td- We also demonstrate how the derived intrinsic 
pulse shapes vary with different PBFs. 

3.1. Simulated Data 

Pulse profiles are generated as single or multi-component 
Gaussians. Scattered profiles are then constructed for a 
given choice of PBF, and noise is added to the scattered 
profile to account for finite signal-to-noise ratios. The 
simulated scattered data are then input to our CLEAN 
algorithm to recover the intrinsic pulse shape as well as 
the best fit model PBF, and the results are compared 
with the original data and PBF used to generate the data. 
Figure 3 illustrates the case where the intrinsic pulse is 
a simple Gaussian and the ISM is described by a simple 
thin slab model. As evident from the figure, the intrinsic 
pulse obtained by application of the CLEAN algorithm is 
in very good agreement with the original data. A signifi- 
cant deviation from the best fit model PBF (a one-sided 
exponential with a characteristic timescale of 40 ms) leads 
to distortion of the resultant deconvolved pulse profile, 
leading to residual scattering (for PBFs narrower than the 
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Fig. 2. — Plots of the various figures of merit used in the CLEAN 
method as a function of the model parameter, tj,. Left-hand side: 
simulated data from Figure 1. Right-hand side: real data from PSR 
J1852+0031 at 1475 MHz also in Figure 1. The vertical dashed lines 
show the true value of for the simulated case and the best value 
for the real data, (a) Nf, the number of data points within 3-<r 
level of the noise, as described in the text; (b) (T ff c /<r ff , the ratio 
of RMS residual noise to the off-pulse RMS noise; (c) the skewness 
parameter, T; (d) the positivity parameter, f r ; and (e) the combined 
skewness + positivity parameter, f c = (7 + f r )/2 

optimal case) or negative going features at trailing edge 
of the pulse (for PBFs wider than optimal case). This 
example thus demonstrates the ability of the technique to 
well characterize the PBF and recover the intrinsic pulse. 

Fig 4 shows another example where the intrinsic pulse 
is modeled as a sum of Gaussians, to mimic a pulse profile 
that is composed of multiple components. As in the pre- 
vious example, the ISM model assumed here is a thin slab 
of scattering material, which can be described by an expo- 
nential PBF with a characteristic time scale Td — 60 ms. 
The recovered pulse is in good agreement with the original 
data. 

3.2. Data from Observations 

Having demonstrated the feasibility of the technique us- 
ing simulated data, we now discuss its application to some 
real data. The data used in this paper arc from recent 
observations at Arecibo, where a large sample of high- 
DM pulsars has been studied at multiple frequency bands 
(Bhat ct al. 2002). The sample includes several new dis- 
coveries from the Parkes multibeam survey (Manchester et 
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Fig. 3. — Example illustrating the CLEAN-based deconvolu- 
tion for simulated data and trial PBFs having a range of pulse- 
broadening times. Top frame: the restored and true pulse shapes 
are shown as solid and dotted curves, respectively, along with the 
true PBF (a one-sided exponential with r = 40 ms) used to generate 
the scattered pulse (dashed line). Lower panel: PBFs and resultant 
deconvolved pulse shapes for a range of models that have the correct 
PBF shape but a range of pulse-broadening times, = 30, 35, 40, 45 
and 50 ms, from top to bottom trace. PBFs that are narrower than 
the true PBF result in pulse shapes with residual scattering (r d =30 
and 35 ms), while wider PBFs over-subtract the data, causing neg- 
ative features near the trailing edge of the pulse (t c ;=45 and 50 ms). 

al. 2001) that are visible at Arecibo. Details of the observa- 
tions and results will be presented in a forthcoming paper. 
Figure 5 shows the 1475 MHz profile of PSR J1852+0031 
(B1849+00), whose DM= 680 pc cm~ 3 implies a distance 
D= 8.4 kpc using the model of Taylor & Cordes (1993). 
We attempted deconvolving these data with two widely 
different types of PBFs: (a) a sharp-edged PBF (eq. 3) 
that characterizes a thin slab model, (b) a rounded PBF 
(eq. 7) for a uniform distribution of scattering material. As 
the figure shows, the results are strikingly different; the 
final CLEAN pulse appears approximately triple for an 
exponential PBF (with a characteristic time r c i=225±14 
ms), whereas it is a distinct double for the case of a uni- 
form medium (t<;=121±6 ms). At first glance, either of 
the models seem equally good. However, given that two- 
component profiles are more commonly observed at and 
near this radio frequency (which implies observations are 
more sensitive to a conal geometry than a core and cone 
case), we believe that a double profile is more likely. The 
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Fig. 4. — Simulation example illustrating the CLEAN-bascd tech- 
nique for recovering an intrinsic pulse shape consisting of three com- 
ponents. The CLEANed pulse (solid line) is overlayed on the origi- 
nal pulse (dotted line); the best fit PBF model (exponential with a 
time constant of 60 ms) is shown on the left (dashed curve). 

double profile obtained with a rounded PBF is constituted 
from only 7 CCs, compared to 21 required to fit the data 
with an exponential PBF. As alluded to before, this econ- 
omy in the number of CCs lends weight to the idea that 
the intrinsic profile is double. Very similar results are ob- 
tained for the same pulsar at 1175 MHz (Figure 6a), which 
implies that although the folded profile at 1475 MHz is an 
average from only ~120 pulse periods, the result is not 
likely to be affected by inadequate profile stability. 

More examples are shown in Figure 6, where the data 
span a wide range of signal to noise ratio and degree of 
scattering. The CLEAN method yields satisfactory and 
robust results even for cases where the scattering tail 
is not very prominent, as illustrated by the scattering 
tail of PSR J1902+0556, which extends only to a few 
percent of the pulse period. The 1175 MHz profile of 
PSR J1852+0031 illustrates the other extreme, where the 
scattering tail extends out to nearly half the pulse period. 
For the 430 MHz data on PSR J1901+0331, three differ- 
ent types of PBFs (eqs. 3, 7 and 7) were used for decon- 
volution. As evident from Figure 6c, the observations are 
better accounted for by a thin slab case with exponential 
PBF of T(j=60±3 ms or by a uniform medium with a PBF 
of characteristic time t c (=31±2 ms, rather than an inter- 
mediate case such as a thick slab (eq. 7). The residual 
asymmetry in the deconvolved pulse is most likely intrin- 
sic to the pulsar beam, though it is possible that the true 
PBF for this line of sight differs from what we have con- 
sidered. Detailed interpretations are deferred to a forth- 
coming paper where we will consider multifrequency data 
and additional forms for the PBF. 

3.3. Estimates of Pulse-broadening Time 

The CLEAN-based approach we outlined in this paper 
provides a general tool that can be applied to a wide vari- 
ety of measured pulse shapes and degrees of scattering. For 
pulse shapes that show significant scattering tails, applica- 
tion of the frequency extrapolation approach as outlined in 
§2.1 can be expected to yield satisfactory results. We have 
applied the CLEAN and frequency-extrapolation methods 
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Table 1 

Pulse-broadening Times: Frequency Extrapolation vs CLEAN a 



F-Extrapolation CLEAN 



PSR 


DM 


Distance* 3 


Frequency 


' a 




T,l 






(pc cm~ 3 ) 


(kpc) 


(MHz) 


(ms) 


(ms) 


(ms) 


(ms) 


J1852+0031 


680 


8.4 


1175 


367 


75 


487 


73 


J1852+0031 


680 


8.4 


1475 


190 


30 


225 c 


14 


J1853+0546 


197.2 


4.9 


1175 


12 


2 


13 


2 


J1853+0546 


197.2 


4.9 


1475 


6 


1 


6 


1 


J1855+0422 


438.6 


10.1 


1175 


25 


15 


24 


4 


J1901+0331 


401.2 


7.8 


430 


57 


3 


61 


2 


J1902+0556 


179.7 


3.9 


430 


15 


3 


14 


1 



a Fits, except where noted, used a pulse broadening function (PBF) that is a one-sided 
exponential, appropriate for a thin screen having a square-law structure function. 

b Derived from DM and the Taylor & Cordes (1993) model for n e in the Galaxy. 

c Data are equally consistent with a rounded PBF (e.g. for a uniform medium) yielding 
Td = 121 ± 6 ms and a double pulse (Figure 5) instead of a merged triple pulse. 



to strongly scattered pulses, such as those in Figure 6, and 
the results are tabulated in Table 1 . The data span a wide 
range of DM, signal to noise ratio and degree of scattering. 
In the frequency extrapolation method, the uncertainty in 
Td is estimated as the model parameter range where x 2 is 
unity above the minimum. For the CLEAN method, a sim- 
ilar estimate can be defined based on the figure of merit 
parameter f r , which is essentially an indicator of the av- 
erage deviation of the residuals from the noise level (nor- 
malized to the noise variance). For instance, if a model 
results in a value for f r that is unity above that for the 
best fit PBF, then the pulse is over-subtracted by a s (on 
average). Note that the best fit model is characterized by 
a minimum in the combined figure of merit parameter f c . 

The CLEAN-based pulse-broadening measurements in 
Table 1 for PSRs J1852+0031 and J1853+0546 indicate 
scaling laws Td cx v~ x with x = 3.45±0.7 and x = 3.4±1.0, 
respectively. These scalings may be compared with those 
expected from a Kolmogorov wavenumber spectrum with 
a negligibly small inner scale [x = 4.4) or from a medium 
with a square-law structure function (x = 4) . The empiri- 
cal intervals for x include both cases, though the centroids 
of the intervals are certainly biased below the expected 
values. We defer a detailed analysis to a later paper, but 
we point out here that the inferred scaling index for a 
given pulsar will depend on the PBF used and that the 
PBF shape may also be a function of frequency. Ideally, 
the data and the CLEAN procedure will indicate the best 
PBF through investigation of the figures of merit discussed 
above. Both the PBF and the scaling index of Td are useful 
information for understanding the nature and distribution 
of scattering along the line of sight (Cordes & Lazio 2001; 
Lohmer et al. 2001; Boldyrev & Gwinn 2002). For mul- 
tifrequency data like these, the CLEAN-based approach 
we have outlined here can be extended to do a global fit 
(by assuming the same form for the PBF, and Td scaling 
with the frequency as a power law but with an unknown 
index) to yield the pulse-broadening times and the index 



of the power law, along with the intrinisic pulse shapes. A 
detailed treatment of this approach is deferred to a forth- 
coming paper where we report on multifrcquency pulse- 
broadening measurements for a larger sample of distant 
pulsars. 

4. CONCLUSIONS 

We have developed and applied an algorithm for de- 
convolving interstellar pulse broadening from radio pulses. 
The technique realizes the two-fold objective of recover- 
ing the intrinsic pulse shape and determining the pulse- 
broadening function that describes the scattering process. 
The method parallels the CLEAN algorithm used in syn- 
thesis imaging applications, but differs from CLEAN in 
that the best pulse broadening function is not known in 
advance. Additionally, we perform a true deconvolution, 
where the resolution of the deconvolved pulse is limited by 
the instrumental resolution rather than the pulse broad- 
ening function. 

Application of the technique is demonstrated through 
use of simulated data as well as data from Arecibo obser- 
vations. Unlike the conventional frequency extrapolation 
approach, the method makes no prior assumptions about 
the intrinsic pulse shape, and can therefore be applied to 
a wider variety of pulse shapes and degrees of scattering. 

The derived intrinsic pulse shape depends on the form 
and parameters of the PBF used for deconvolving, and 
further considerations have to be employed to select the 
best fit model. We employ figures of merit that quan- 
tify positivity and asymmetry of the CLEANed pulse, and 
the best fit PBF is taken to be the one that yields mini- 
mal residual scattering and skewness of the intrinsic pulse. 
As in any other nonlinear method, the uniqueness of the 
deconvolved result is not always assured: however, further 
considerations (such as data at different observing frequen- 
cies) can help resolve any such degeneracy. Intrinsic pulse 
shapes obtained by application of this method are likely 
to be closer to reality, and hence should prove to be useful 
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Fig. 5.— CLEANed pulse shapes for PSR J1852+0031 at 1475 
MHz obtained by deconvolving with two different models for the 
PBF. The results are strikingly different. Upper panel: a one- 
sided exponential PBF, corresponding to a thin-screen model with 
a square-law structure function, yields t ( j=225±14 ms and an ap- 
proximately triple pulse with merged components. Lower panel: 
a more rounded PBF, corresponding to a uniform medium with a 
square-law structure function, yields 7^=12^6 ms and a distinct 
double pulse shape. The differences result from the fact that in the 
top panel, the sharp PBF requires broad structure in the decon- 
volved pulse in order to match the measured pulse. For the bottom 
panel, the PBF is more rounded and the deconvolved pulse can have 
sharper structure and still be consistent with the measured pulse. 



for the study of pulsar emission properties. A paper is in 
preparation which applies our method to a large body of 
multi- frequency pulsar observations. 
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Fig. 6.— The measured and CLEANed pulse profiles for four pulsars: (a) J1852+0031 at 1175 MHz, (b) J1853+0546 at 1475 MHz, (c) PSR 
J1901+0331 at 430 MHz; The measured pulse can be accounted for by either a thin slab model with a PBF of characteristic time, r d =60±3 
ms (dashed curves), or by a uniform medium with a PBF with Td=31±2 ms (solid curves), and (d) PSR J1902+0556 at 430 MHz. The best fit 
ISM models (PBFs) are shown as dashed curves (panels a, b and d), and the CLEANed pulses are shown as thick solid curves. The negative 
dip preceding the pulse in (a) is an artifact that appears for several choices of PBF and tj. 



